Goos-Hanchen and Imbert-Fedorov shifts for bounded wave packets of light 
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We present precise expressions of the spatial and angular Goos-Hanchen and Imbert-Fedorov 
shifts experienced by a longitudinally and transversally limited beam of light (wave packet) upon 
reflection from a dielectric interface, as opposed to the well-known case of a monochromatic beam 
which is bounded in transverse directions but infinitely extended along the direction of propagation. 
This is done under the assumption that the detector time is longer than the temporal length of the 
wave packet (wave packet regime). Our results will be applied to the case of a Gaussian wave packet 
and show that, at the leading order in the Taylor expansion of reflected-field amplitudes, the results 
are the same of the monochromatic case. 



<Z3 ■ I. INTRODUCTION 

O 
• i-H 

The problem of the reflection of an electromagnetic field upon a planar surface is an old and well-known problem 
Q , of optics which is usually treated by means of the so-called Fresnel reflection (and transmission) coefficients. These 
are typically given for the simplest case of impinging plane waves solely [IHi. For beam-like field configurations this 
is not a problem since any electromagnetic field of arbitrary shape can be decomposed into plane waves. Thus, the 
reflected field can be evaluated by first reflecting each of its plane wave components individually and then sum them 
£>-j" back. However, non-specular reflection effects mayoccur and be observed [5|-|7[. Among all these effects, the most 
Jin . famous ones are the so-called Goos-Hanchen (GH) |8l-fl2j and Imbert-Fedorov (IF) fl3l - [l6j shifts, the former occurring 
in the plane of incidence of light while the latter occurs in the plane orthogonal to the plane of incidence. Detailed 
theoretical predictions and experimental observation of both GH and IF shifts have been given in the past decades for 
— ■ different types of surfaces (dielectric, metallic [13, EH, multilayered [HlHil, et cetera) and for different beam shapes 



t> | (Hermite-Gauss [2JJ, l22j, Laguerre- Gauss [23j|, Bessel [2J| et cetera). Despite the vast amount of literature on the 



subject, however, there is almost no clear and simple theory for the non-monochromatic case, namely for beams of 
light which have not only finite transverse size, but also finite temporal extension (wave packets). A notable exception 
is represented by the work of Gragg (?] whose treatment is however limited to the case of scalar waves. 

In this work, starting from the work done in Ref. [3] and following the formalism of beam shifts as displacements 
of the center of mass of the energy density distribution of the electromagnetic field upon reflection introduced in Ref. 

[U 

l26j . we extend the concepts of Goos-Hanchen and Imbert-Fedorov shifts to the non-monochromatic case, firstly 
by generalizing the usual monochromatic formulas to the non-monochromatic case, and then by applying these new 
results to the simple case of a Gaussian quasi monochromatic wave packet. 

This work is organized as follows: in Sect. 2 the notation used throughout the manuscript is established and the 
plane- wave expansion of a non-monochromatic field is given. This will be used in Sect. 3 to calculate the center of mass 
of the field energy density distribution, which will allow us to calculate the shifts. Then, in Sect. 4 a brief recall on 
how to calculate the reflected field amplitudes is given, and the general expression for the non monochromatic center 
of mass is derived explicitly. In Sect. 5, the application of this formalism to the simple case of a quasi monochromatic 
Gaussian wave packet is discussed. Finally, in Sect. 6 conclusions are drawn. 
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II. FOURIER REPRESENTATION OF A NON-MONOCHROMATIC FIELD 



A. Notation 

Throughout the manuscript, the following conventions are understood. According to Ref. [13], we indicate with 
capital letters the transversal components of a three dimensional vector, so that we will write r = {x, y, z} — {R, z}, 
being R = {x,y} for the position vector and k = {k x , k y , k z } = {K,fc z }, being K = {k x ,k y } for the /c-vector. 
Moreover, we define r = |r|, R = |R|, k = |k| and K = |K|. We also define the gradient operator in the A:-space as 
follows: 

-[— — —) = !-— 1 

9k \ dk x ' dk y dk z j \dK dk z / v ; 

Finally we assume the following condition on the summation indexes: Greek indexes can assume only two values, 
namely /i, v 1 r, • • • £ {1, 2} while Latin indexes assume three values: {1, 2, 3}. 



B. Local basis 



Let us consider a Cartesian reference frame (Oxyz) associated to the orthonormal basis {x, y, z}. Let ko = fcoz be 
the central wave vector of a wave packet propagating along the z-axis, and k = fc x x + k y y + k z i denotes the generic 



wave vector with k — yk^ + ky + k z k . Given k = k/fc and the unit vector z parallel to the z-axis, it is possible 
to build a complete orthogonal basis e.j (k) out of the standard Cartesian basis by means of a rotation around the axis 

' (2) 



|z x k| 



by the angle 9 — arcsin(if/fc) between k and z [28j], where "x" denotes the ordinary cross-product in R 3 . By 
representing such a rotation with the Rodrigues' formula [29[ 

R(9,n) = e eM = / + Msin0 + M 2 (l - cos (9), (3) 

where / is the (three dimensional) identity matrix and the matrix elements of M are given by: 

3 

M v = -^2 £ *jk(n)k, (4) 
fc=i 

where e^-fc is the completely antisymmetric Levi-Civita symbol. By doing this operation the local basis e^k) can be 
explicitly written as ed(k) = R(6,h)5t, 62 (k) = R(8,h)y and §3 (k) = R(6,h)z = k. It is worth noticing that this 
basis constitutes a local basis attached to the field itself, and that 63 (k) points in the direction of the wave vector k, 
while the other two basis vectors ei^k) are orthogonal to k and they are natural candidates for describing the field 
polarization These basis vectors have the following properties: 

e M (k) • e„(k) = 5^, (5a) 
ei(k) x ej(k) = e^e^k), (5b) 

( k x e M (k)) • f k x ^(k)^ = 5^, (5c) 

where the dot symbol indicates the ordinary scalar product in R 3 and summation over repeated indexes is understood. 
Note that with a suitable unitary transformation [l| it is possible to relate this basis with other basis commonly used 
in optics, like the helicity basis. 
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C. 3D Fourier representation of the electromagnetic field 

Having established the form of the local basis attached to the (incident) beam field, we can now write the electric 
and magnetic fields in terms of plane waves as follows [l| : 

2 

E(r,t) = / ^*A |i (k)e M (k)e i ( k - r - t,rt >+c.c J (6a) 

2 

cB(r,t) =^2 d 3 kA^{k)[k x e^k)]^^-^ +c.c, (6b) 

where = c|k| = c^J k% + + k%, d 3 k = dk x dk y dk z and "c.c." stands for complex conjugate. The complex amplitude 

An(k) contains informations on both the envelope and on the polarization of the single spectral component of the 
field, and can therefore be written as A M (k) = A(k, ko)a M (k), where A(k, ko) is the scalar distribution of the plane 
waves components of the wave packet around the central wave vector ko and a^(k) contains informations about the 
field polarization. 

As a final remark to this section, it is worth noticing that in the integrals that appear in Eqs. ([6]), k x , k y and 
k z are independent variables. In this general formulation, no hypothesis on the nature of the field (whenever the 
field is homogeneous or contains evanescent contributions) is made, and consequently no relation between the various 
components of k has been established yet. 



III. CENTER OF MASS OF THE FIELD ENERGY DENSITY 

In this section we will generalize the formalism introduced in Ref. [26] for the calculation of monochromatic beam 
shifts. The essential ingredient is the center of mass of the energy density distribution, defined as follows: 

[ dt [ d 2 RRH(R,z,t) 

(R) = T T > (7) 

J dt j d 2 RH(H,z,t) 

where % = £o(E 2 + c 2 B 2 )/2 is the electromagnetic field energy density and T is the integration time of the detector. 
Note that, compared with the expression for a monochromatic field (33|, here the presence of the integral on time 
is imposed by the non monochromatic nature of the field. In order for the detector to be able to detect the whole 
energy of the wave packet (WP), the integration time T must be chosen in such a way to be larger than the length 
in time of the WP. If we call £ the width of the WP along the z direction, with (ko <C 1 , T must be compared 
with £/c. Thus, if the integration time is much longer than the temporal length of the WP (i.e., cT ^> £), we are 
permitted to let T — » oo in the integrals in Eqs. ([7]). Conversely, if the integration time is comparable with the WP 
timescale (namely if cT < £) , we are not allowed to extend the integration extrema of the time integral to infinity, 
as the detector can detect only a fraction of the energy of the WP. In this case, then, the role of the detector will be 
actively important, since the choice of the integration time can drastically modify the properties of the detected field. 
In this work we will focus our attention on the first case (the so-called wave packet regime), leaving the other case 
(the beam regime) to future analysis. 

With this in mind let us derive the expression for the center of mass (R) . The electromagnetic field energy density 
is given by |2j 

U = U e + U b = j [E(r, t) • E(r, t) + c 2 B(r, t) • B(r, t)] . (8) 

By substituting the expressions of the electric and magnetic fields in terms of their spectral component as defined in 
Eqs. ^ we obtain, after a long but straightforward calculation, the following expression for the center of mass: 

+ z^\e u (k)A u (k) 
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E/dffc^ck)fi;(k).{-i 

X) / d 3 fc£iA M (k)r 



<R> = ; — , (9) 
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where k/k z > plays the role of a weight function that approaches the value 1 for well collimated beams such 
that (k 2 . + k 2 )/k 2 <C 1. In obtaining this expression, some care has to be put when integrating over the transverse 
coordinates and over time. While the integration over d 2 R produces a Dirac delta on the transverse fc-vector K, 
leaving an explicit dependence on k z , namely 

d 2 R e *( k - k ')- r = (2tt) V^- fe ^M(K' - K), (10) 

the time integration has to be analyzed more carefully. To this aim, we will now show how to pass from Eq. ([7]) to 
Eq. ©. To do that, we will concentrate on the contribution of the electric part of the electromagnetic energy density 
to the center of mass ([7]), as the magnetic contribution can be obtained using similar calculations. Let us start by 
evaluating the denominator of Eq. ([7]). By substituting the expression of the electric field given by Eq. (|6a| into the 
expression of the electric part of the electromagnetic energy density that appears in the denominator of Eq. ^ we 
have: 

I, = I dt 

£o_ 

A./; 

X 



J J d 2 RH e (TL,z,t) 
2(2 ^)3 E / d " k J d 3 k'A;(k)A x (k')e;(k)-e x (k') 

Jdte^Jd 2 R^ + , C , (11) 

where the "fast-oscillating" terms leading to 8{uj + lu') — > have been omitted. This expression can be evaluated by 
using Eq. (fTUf and the standard properties of the Dirac delta function [32j yielding to: 

J dfe-*(«-«')t = 27 n5(w' - u) = (^j ^- [s(k' z - k z ) + S(k' z + k z )] , (12) 

where the two Dirac delta functions 5{k' z ± k z ) correspond respectively to the forward and the backward propagating 
solution. Without any loss of generality, we assume to consider only forward propagating waves (for that we can write 
\k z \ — k z ). Substituting these expressions into the previous equation leads to the following result: 

J « = £ -fY,f d3k J d 3 k>^A;(k)A x (k')e;(k)-e x (k>) 

2 r k 

x eW'-WstK' - K)S(k' z -k z ) = ^J2 d3k T I^WI 2 ' ( 13 ) 



c 



For the denominator, then, the presence of the time integration is not critical, since the only effect is to introduce a 
Dirac delta in the field frequencies 8{u — u) that is easily transformed into a Dirac delta in the k z component of the 
fc-vector times the factor k/k z . 

For the numerator, the calculations need a bit more of attention. First of all, for the sake of simplicity, let us 
consider only the contribution given by xE(r, t), since the contribution of the term ?/E(r, t) together with its magnetic 

counterpart can be straightforwardly calculated with similar techniques. By calling <?R x H e (R, z, t) 

and by remembering that (when integrating in d 3 k) the presence of the multiplicative term x can be turned into a 
derivative with respect to k x by means of part integration we obtain the following result: 



so 



2 ^ 3 tr 



J d 3 k J d 3 k' J dt J J d 2 i?A;(k')e;(k')e i ( k ' r - w ' t ) 




A x {k)e^k)e l ^ r -^ 



+ c.c\, (14) 



where, again, the "fast oscillating" terms have been omitted as they give zero contribution. Note that the calculation 
of the derivative brings four terms. While the first three of them (the derivative of the amplitude, the derivative 
of the local basis and the derivative of the spatial plane wave respectively) have a time dependence of the form 
exp [— i(u) — oj)t] that leads, once integrated, to a Dirac delta in the field frequencies 6(ui — w), the fourth term contains 
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the derivative of co with respect to k x (since u = c|k|). This brings down a factor of t from the time exponential that 
gives rise (once integrated in time) to a derivative of the Dirac delta with respect to u, i.e., dS(uj' — u))/du). 

Since u> is not an integration variable but it depends on |k| via the vacuum dispersion relation lu = c« / k x 



kl 



ck (and equivalently ui' = c^J k% + k^ + k' z 2 = ck'), it is more convenient to express this Dirac delta as a function of 

the variable k z , in order to relate it to one of the integration variables and use it lately to eliminate the integration 
in k' z . By using Eq. (|T0|) together with the following formula [32j 



dtte 



kk' 



_d_ 



_d_ 

dk z 



kk' z + k'k 

it is possible, after some straightforward algebra, to obtain the following result: 



5(u'-u), 



(15) 



lex — 



^E J >^W<w 



A,// 

el(k) 




i\dk x k z dk z W 



(16) 



The same result can be obtained for the y-component by making the substitution k x —¥ k y . 

In order to make Eq. @ more easy to understand and manage, it is possible to exploit the quantum notation for 
wave packet and beams of light introduced in Ref. [34j to rewrite Eq. in a more appealing way by separating the 
contributions to the center of mass of the energy distribution in three parts as follows: 



(R) 



(R)s + (R) 



z(K) 



(17) 



where the first term accounts for spatial shifts, the second term for the contribution of the Berry connection [34| and 
the third z-dependent term accounts for the angular shift. Their explicit form is given by: 



1 


' d 


K 


( d 1 


i 


dK 


k z 


V dk z k z 



\A L 



(R 



^E^IirV(k)l^), 

R )^ = ^E(^ 



k 



K 



where 



A^(k) = e;(k) 
is the Berry connection expressed in the basis e^(k) and 



1 / d K d 

i \ 9K k z dk z 



e,(k) 



^=E^i r i^ 



(18a) 
(18b) 
(18c) 

(19) 
(20) 



is the normalization factor (i.e., the denominator of Eq. ((9|)). This is the first main result of this work. As can be 
seen from the previous equations, a first principal difference with respect to the monochromatic case is the presence 
of the weighting factor k/k z . This term, that comes from the conversion of the Dirac delta 8(lj — lu) (yield by time 
integration) into a Dirac delta S(k' z — k z ), takes into account the fact that the component k z of the wave vector k 
can now independently take any value and it is not anymore constrained to be dependent on k x and k y as in the 
monochromatic case. 

The other important difference with respect to the monochromatic case is the presence of the last factor (the term 
l/k z in the round parenthesis of Eq. (|18ap ). In order to understand the origin of this term, let us consider explicitly 
the calculation of Eq. (I18a[) . 
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Starting from Eq. (|T4|) . after having performed the integration in the transverse spatial coordinates and the two 
time integrations (one that gives 5(u' — u>) and the other one that gives dS(uj' — to)/du>') we obtain the following 
expression (we drop the constant factors and the summation over \x since they are of no relevance for this discussion): 



d A k- 



•V k >( -i^? ) + MV [ li? 



1 d \ 
~ldKJ 




(21) 



where the z-dependent part and the term that contains the Berry connection have been ignored since the last term 
only originates from these two terms. By noting that (since |A M (k)| 2 is real) 



where kj £E {k x , k y , k z }, the previous term can be rewritten as follows 



1 / d K d 
7 V <9K ~ ~kZ~dk 2 



d 6 k- 



11 d K d 
~i\ dK ~ T z ~dk z 



(22) 



(23) 



The first integral gives the usual term that depends on the linear momentum, while the second term, after a partial 
integration, accounts exactly for the extra term in Eq. (|18al) . This extra term, then, is due to the non monochromatic 
nature of the field, since it involves, ultimately, the derivative of the weighting factor k/\k z \, namely 



I d K d 

i \ dK. fe, dk. 



(24) 



In the monochromatic case, in fact, since the weighting factor is constant, this term vanishes. 

IV. CALCULATION OF THE REFLECTED AMPLITUDES 

The geometry of the problem is depicted in Fig. [1] We indicate with k the representation of the A;-vector in the 
reference frame attached with the incident field, while with k we intend the representation of the k- vector in the frame 
attached to the reflected beam [7|. Upon reflection, the relation between these two A:- vectors is the following: 



k = k x x + k y y - k z z. 



(25) 



The following rules for scalar product also trivially hold: k • r = k • f and k • r = k • f . Having written the electric and 
magnetic fields in terms of their 3D Fourier spectral components, all that is needed for calculating the amplitudes of 
the reflected field is to know how the single plane wave components reflect. Following Ref. [|[ it is easy to show that 
in order to calculate the reflected field amplitude it is enough to do the following substitution: 



A M (k)e M (k)e 



ikr 



A,(k)r,(k)e,(k)e^ = A'(k)e^(k)e 



_ik-r 



(26) 



where r M (k) is the Fresnel reflection coefficient for the p-polarized waves (or TM polarization, that is plane waves 
with the electric field parallel to the plane of incidence) for fj, = 1 and for s-polarized waves (or TE, that is plane 
waves with the electric field orthogonal to the plane of incidence) for fi — 2. The explicit expression for the Fresnel 
reflection coefficient is given by 



r P (k) 



ek z 



(27a) 



where e is the (generally complex- valued) relative dielectric constant ei /e2 of the two media that constitute the interface 
and fc* = Jek—k^. — k% is the ^-component of the wave vector inside the medium. Note that these expressions are 
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FIG. 1: Geometry of the problem. {x;,yi,Zi} is the reference frame attached to the incident wave packet, {x r , yy, z r } is the 
reference frame attached to the reflected wave packet and {x, y, z} is the laboratory frame. 9 is the angle of incidence. 



conceptually different from the one normally given in textbooks l3( for the monochromatic case. Here, the term 
k is no longer constant and it cannot be factorized outside the square root anymore. Thus, in the expression of the 
reflection coefficients for a non monochromatic fields, an intrinsic dependence of the coefficient on the frequency uj (or 
analogously on the k z component of the wave vector) is present. 

The expression for the center of mass of the energy distribution in the reflected frame is then equal to the one 
obtained in the incident frame (Eq. (|9])) with the only difference that now we have to consider as weighting functions 
the amplitudes of the field in the reflected frame, namely we need to make the substitution \A^) — > {A'^) and 
(A^l —> (A'\. Then, following the prescription given in Ref. (26|, given a center of mass distribution (R)(z) like the 
one calculated in Eq. ([9]) it is possible to compute both the spatial and angular Goos-Hanchen and Imbert-Fedorov 
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shifts by means of the following relation: 



A sp 



^IF 



(y) 



z=0 



z=0 



A ang 
^GH 



A ang 
IF 



d_ 

dz 
d_ 

dz 



(y)- 



(28a) 
(28b) 



Note that in our case, by using the expression given by Eq. (1171) for the center of mass distribution and this definition 
of the two shifts, the expression of the center of mass can be put into a direct relation with the Goos-Hanchen and 
Imbert-Fedorov shifts thanks to the following relation: 

(R) = x Aqh (z) + y Aif (z) , (29) 

where, according to Eqs. ([2"81 

AghW = (x) s + (x) B + z(x) A = A^ P H + zA^f , (30a) 

A lF (z) = (y) s + (y) B + z{y) A = + zA^ s . (30b) 

Note that this expression that connects the center of mass with the shifts is very general (namely, it is beam-shape 
independent). The explicit expression of both Goos-Hanchen and Imbert-Fedorov shifts can be calculated once the 
explicit shape of the beam is known. 



V. APPLICATION TO A GAUSSIAN WAVEPACKET 



In this section we will now apply the formula given above for the case of Gaussian wave packet. We assume 
that it is possible to write this wave packet as a narrow distribution of A:-vectors around a central value ko (quasi 
monochromatic approximation). We also define in the reference frame {OxiyiZi) attached to the incident beam the 



vector f = a p Xi + a s exp (i(f>)yi that accounts for the beam polarization, with a p ,a s ,< 
expression of the polarization function a M (k) is then given by 

,1-1/2 



€ R and al + 



1 - 



If • kl- 



<(k)-f, 



and the spectral amplitude ^4(q) is assumed to be a Gussian distribution, namely 



A(q) 



(detV) 
^74" 



1/4 



exp 



-iq ■ r - -q • (Vq) 



1. The 



(31) 



(32) 



where we defined q = k — ko and ro accounts for some possible initial transversal beam displacement. The matrix V 
is given, in the frame (OxiyiZi), by 



l/o-i 
V = I I/0-2 
I/03 



(33) 



being <ri 2 the width of the Gaussian wave packet in the transversal plane {k Xil k yi } and 03 the width of the Gaussian 
wave packet in the longitudinal direction k Zi . 

By making the change of variables q = k — ko in Eqs. (|18p and by expanding in a Taylor series all the terms that 
appear in the integrals as a function of ko+q (namely the reflection coefficients r^ko+q) and the polarization functions 
a M (ko + q)) up to the first order for the spatial and the Berry shift and up to second order for the z-dependent part, 
and by representing the (generally complex- valued) reflection coefficients using the polar representation of a complex 
number as 



1 (k) = |r M (k)|e^( k ) = i? A1 (k) £ 



*<M k ) 



we obtain the following results 



E 



"dk x 

1/2 

y (w p w s ) 1 cot 9 si 



A w p a 2 s + w s a 2 . 
y cot sin < 



E 



X Sfw^ 



sm( 
01ni? p 

dk x 



y s 2 w^ 



d\nR„ 



Ok 



v 



w s a 



cot 6 cos < 



CLpd s 



(34) 

(35a) 

(35b) 
(35c) 
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where we have defined the weighting factors (namely the fraction of energy contained in each polarization) as 



«2 r?2 



• alRl + a^RV 



(36) 



and the factors Sf (i — {1, 2, 3}) are defined in the following manner: 



Si 



y/2' 



(37) 



where 9i is the angular divergence of the wave packet along the direction Xi (with x\ = x, 22 = y and X3 = z). 

Summing all these contributions we obtain the final expression for the center of mass of the energy density distri- 
bution, that can be written as 



(R)= £ [xA™(z)+yA^(z)] 



(38) 



where the Goos-Hanchen and Imbert-Fedorov shifts are given by the following formulas: 

'dtp ^ 



KM 



b2 dk„ 



w p a 2 s 4- w s a 



„ . w p a 2 s - w s a 2 

— sm (p — z cos < 

dpd s 



(w p w s ) 1/2 sin(0 -<t> p -4>s) 



cot< 



(39a) 



(39b) 



Note that in writing these formulas we retained the terms depending on the derivative with respect to k y , often 
neglected because they are zero for Fresnel coefficients of a planar interface. Here, however, we kept them for the sake 
of completeness. 

Note, moreover, that the results presented in this section give the same expression of the GH and IF shift of a 
monochromatic beam of light, despite they refer to the non monochromatic case. This fact is due to the fact that we 
stopped our analysis to the first order (in Taylor expansion) for the spatial and the Berry term and second order for 
the angular term, while, as stated in Ref. non monochromatic corrections to the GH and IF shifts should manifest 
(as an explicit dependence of the formulas above on the longitudinal width S3) starting from the fourth order in the 
expansion with respect to k — ko. 



VI. CONCLUSIONS 



In conclusion, in this work we extended the concepts of Goos-Hanchen and Imbert-Fedorov shifts to the case of a 
non-monochromatic beam of light, namely a bounded wave packet. We shown that the expression of the center of 
mass of the energy density distribution can be easily written in a form similar to the monochromatic case (where 
the presence of the weighting factor k/\k z \ is the signature of non monochromaticity) and we applied these results to 
the simple case of a quasi monochromatic Gaussian beam, showing that the non monochromatic nature of the beam 
does not play a role up to the first and second order terms of a suitable perturbative expansion, for the spatial and 
the angular part, respectively. These results may appear simple at a first glance, but in reality they are not, since 
conceptually it is not obvious that at the lowest order the non monochromatic nature of a field does not play any role 
in determining its shifts upon reflection. Moreover, since the expression of the center of mass for a non monochromatic 
wave packet of light (Eq. ([9])) is sensibly different from its monochromatic counterpart, it is not trivial at all to expect 
that these modification do not count or exactly compensate. Moreover, the presence of the weighting factor k/\k z \ 
makes the integrals not solvable anymore analytically if this factor is not expanded in series. This is a signature of 
the perturbative action of the non monochromaticity on the process of reflection, and it clearly states that the major 
contribution comes from the zeroth order, i.e. from the monochromatic results. In future works, however, we intend to 
study deeper this problem, trying to understand under which conditions (and at what order) the non monochromatic 
corrections start to play a role. 
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